Non-equilibrium dynamics of interacting Fermi systems in quench experiments 
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We describe the dynamics of two component non-interacting ultracold Fermions which are initially in thermal 
equilibrium and undergo a rapid quench to either the repulsive or attractive side of a Feshbach resonance. The 
short time dynamics is dominated by the exponentially growing collective modes. We study the Stoner instability 
and formation of ferromagnetic textures on the repulsive side, and the pairing instability towards BCS or FFLO- 
like states (determined by the population imbalance) on the attractive side. In each case, we evaluate the growth 
rate of unstable modes and predict the typical lengthscale of textures to be formed. 



Non-equilibrium dynamics of strongly interacting quantum 
many body systems is a subject of great interest with diverse 
applications in different fields, from condensed matter physics 
to high energy physics to cosmology. However, compared 
to equilibrium physics, it lacks deeper understanding from 
both theoretical and experimental standpoints. Recent devel- 
opments with ultracold atoms [ 1 ] have opened a new window 
to study strongly interacting model Hamiltonians. Although 
the focus has been mostly on studying the equilibrium phases 
of these models, these systems also present a unique oppor- 
tunity to probe the intrinsic non-equilibrium dynamics in the 
strongly correlated regime [2|. The low energy scales in the 
system (typically ~ kHz) and an almost ideal insulation from 
the environment results in comparatively larger timescales of 
dynamics than in the corresponding condensed matter systems 
and therefore, make it easier to follow the intrinsic quantum 
evolution of the system in real time Q. 

The inherent tunability of interactions in these systems al- 
lows us to study the dynamics of the system after a quench 
(rapid change) of the Hamiltonian parameter. In fact, quench 
dynamics has already been studied in cold atom experiments 
that probe dynamics and pattern formation in spinor BEC Pfl . 
collapse and revival in optical latices Q, and "Bose novae" 
in BEC's with attractive interactions (6). Our primary moti- 
vation for this study comes from recent experiments on cold 
two-component Fermions in which interactions are tuned to 
the strong repulsive regime during the course of the experi- 
ment [7]. It was found that for sufficiently strong coupling, 
the atom loss rate, the size of the atomic cloud and the dis- 
tribution of the kinetic energy of the atoms show signatures 
of a transition to a ferromagnetic state. Originally, these ex- 
periments were interpreted as probing the equilibrium Stoner 
transition (8). On the other hand, the time scales for tuning the 
interactions between small and large values were limited by 
the formation of molecules on the BEC side of the Feshbach 
resonance. In this paper we show that the finite experimental 
time scales are such that we can also analyze some aspects of 
the experiments from the point of view of quench dynamics. 

In this Letter, we focus on the initial dynamics of the col- 
lective modes in Fermionic systems following a rapid quench. 
Within linear response theory, we find that for quenches across 
a critical interaction strength U c , some collective modes be- 



come dynamically unstable. We denote the growth rate cor- 
responding to an unstable mode with wavevector q by A q . In 
general, the lengthscale of the textures formed in the quench 
is bounded by the smallest unstable wavelength q^. Further- 
more, if the most unstable mode (i.e. for which A q assumes 
its maximum) has a finite momentum, q mdx , then the order pa- 
rameter will have large spatial modulations at the correspond- 
ing wavelength q~.} x . We study two typical cases of dynami- 
cal instability in Fermi systems starting with a non-interacting 
gas: (i) quench to the regime of strong repulsive interactions 
which leads to the Stoner instability and ferromagnetic tex- 
tures; (ii) quench to the regime of weak attractive interactions 
and the BCS pairing instability. The order parameter in the 
former case (magnetization) is conserved, while in the latter 
case (wavefunction of the coherent pairs), it is not conserved. 

Our main results are: (a) for the Stoner instability in 3D, the 
wave-vector of the most unstable mode scales as q mdx ~ it 1 / 2 
while its growth rate A max ~ u 3 / 2 , where u = Uf/U c — 1 and 
Uf the final coupling after the quench, (b) In 2D, at T = 0, 
the most unstable mode has growth rate ~ u and occurs at 
\f2kp, where hp is the Fermi wave- vector of the gas. How- 
ever, the system is very sensitive to thermal fluctuations and 
at finite temperatures, g max scales with u 1 / 2 T 1 / 2 exp(eF/2T) 
while the growth rate ~ u q max . This scaling reverts back to its 
T = form once q m . dx ~ \f2kp. (c) In the BCS case (in 3D), 
we find that for small initial spin imbalances, the fastest grow- 
ing instability is towards formation of coherent Cooper pairs 
with center of mass momentum q = 0. Beyond a critical po- 
larization, the fastest growing modes correspond to pairings at 
a finite wavevector, i.e. towards FFLO-like states |9, 10 1. Be- 
yond a 2 nd critical polarization, all modes become stable. We 
note that the change in the character of dynamics does not nec- 
essarily correspond to crossing a thermodynamic phase tran- 
sition during the quench; e.g. although there is no ferromag- 
netic transition in 2D at finite temperatures, quench dynam- 
ics still involves unstable modes which exhibit exponential 
growth. 

Formalism- We consider a system of interacting Fermions 
described by the Hamiltonian: 

H = Jacket,, + U(t) J d d v^cl iCri c rV (1) 
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FIG. 1 : Properties of growing collective modes in the Stoner insta- 
bility in 3D (top row) and 2D (bottom row), (a) and (c): Growth 
rate A q as a function of wavevector q for T = (thick lines) 
and T = 0.1 Tf (dashed lines) and several values of the dimen- 
sionless coupling u = Uf/U C (T) — 1 = 0.01 (bottom pair of 
curves), 0.03, 0.05, . . . , 0.13 (top pair of curves). The insets show 
the dynamics phase diagram as a function of final interaction strength 
Uf /U c {0) and reduced temperature T/Tf- "PM" corresponds to the 
case where all collective are damped, while "FM" corresponds to the 
existence of an unstable collective mode, (b) and (d): The most un- 
stable wavevector g max and the corresponding growth rate A max vs. 
the dimensionless coupling u for T = (solid lines, both 2D and 
3D), T = 0.2 T F (dotted lines, 3D), T = 0.1 T F (dotted lines, 2D). 
The dashed lines in (b) are the approximate analytical solutions. 



where c^c^) are the Fermion creation (annihilation) opera- 
tors with spin a, £ker = k 2 /2m — m is the mass of the 
Fermions, /Lt CT are the chemical potentials, and U(t) is the time 
dependent coupling. U > corresponds to the repulsive case 
(Stoner instability), while U < corresponds to the attractive 
case (BCS pairing instability). We focus on the instantaneous 
quench limit, in which the coupling U changes from a negli- 
gible initial value [7j to a final value Uf at time t = 0. Ut 
is chosen to lie beyond the critical interaction strength U C (T) 
for the transition in question. The instantaneous quench ap- 
proximation is valid if the time of ramp is much shorter than 
the growth rate of the instability (will be discussed later). In 
rapid quenches the distribution of fluctuations present in the 
system does not change during the quench, and is given by 
the equilibrium distribution at the initial value of interaction 
U. We are interested in the dynamics of this initial distribu- 
tion after the quench. This dynamics is dominated by unstable 
collective modes which grow exponentially in time and lead 
to formation of textures. The equation of motion for the col- 
lective coordinate $ q is given by: 



iSe$ q (t) = [$«,(*), fl], 



instability, while for 
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for the Stoner 



instability, <f> q (t) 



Sk c ic+q t-k j.' To probe the dynamics, we couple an ex- 
ternal field Bq(t) to the collective coordinate and look at its 
response to this perturbation (for Stoner, the external field is 
just the magnetic field while for BCS, it is a fictitious pairing 
field). Within the framework of RPA, we obtain: 

= / x2(* - [f(0<MO + £ q (0] , 0) 



where $ q (<) = ($ q (i)) and - t') is the retarded bare 
susceptibility in the initial (unordered) equilibrium state. This 
assumes that the Fermion distribution functions do not change 
during the short time dynamics after a rapid quench. Eq. ([3]) 
may be interpreted as the propagation of the total field (com- 
posed of the internal mean field £/(£')<I> q (i') and the external 
field Bq(t')) forward in time by x q (^ — 0- ^ tne external 
field Bq[t'), as well as $ q (t') vanish for t' < flT], Eq. ^ 
can be explicitly solved for $ q (i): 



*q(*)= / dt'x^t-t'-U^B^t'). 



(4) 



where Xq(^ — t'\ Uf) is the inverse Fourier transform (with 
the integration contour shifted above all poles in the complex 
plane) of 



Xq(w;E//) 



1-UfX^coY 



(5) 



This useful result relates the short time dynamics of the col- 
lective mode after a quench to an RPA-like susceptibility in 
which the bare response functions are evaluated in the ini- 
tial state, but with a coupling corresponding to the final state. 



When Xq( w q; Uf) has poles at i 



iAq in the upper 



half of the complex plane, then fluctuations that occur after the 
quench will grow exponentially in time. In the next sections, 
we make these ideas explicit by studying the specific cases of 
the Stoner and the BCS instability. 

One of the common concerns about the existence of the 
Stoner instability in equilibrium systems is the so called 
Kanamori argument [12|, which says that the effective U 
which enters the Stoner criterion gets reduced by screening. 
We expect that the time scale associated with the osscillations 
of the screening potential is ~ e^ 1 which is much faster than 
A q 1 . We proceed by assuming that there is a phase transition, 
so, the U f that enters Eq.(|5]) should be considered to be the 
renormalized U f . 

Stoner instability- For a quench of the interaction across 
the Stoner instability, the magnetization propagator is given 
by Eq. (|5]), with 



d d k n F (&+ q ) - n F (&) 
(27r) d W -(£k+ q -6c) + *0+ : 



(6) 



^) where X q ( w ) i s m e bare spin susceptibility per spin of a free 
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Fermi gas, and n F is the Fermi distribution function. The 
poles of the corresponding propagator are then found from: 

l-0) X °(fi q + iA q ) = O. (7) 

A line of complex poles with a positive imaginary part A q , 
appears when Uf exceeds U C (T) which corresponds to the 
Stoner instability (see insets of Fig. [T^ and c). A q grows 
linearly for small q, reaches a maximum value, A max at a 
wave-vector g max , and decreases to zero at g cut . Modes with 
q > 9cut are stable. 

We begin by considering the 3D case. At T = 0, the in- 
tegrals in Eq. (|6]i may be evaluated explicitly, yielding the 
familiar Linhardt function, analytically continued to the up- 
per half-plane. We obtain the approximate expression for A q 
at T = 0: A 3D = (2k F /irm)uq - (l/6irmk F )q 3 + 0(q 5 ). 
Fig.[TJi shows A q asafunction of g,forT = 0andT = 0.1 T F 
and several different values of u = Uj/U C (T) — 1. In all 
cases A q=0 = as required by the conservation of total 
magnetization. The values of A max and <? max are plotted as 
a function of u in Fig. [lj>. The analytic expression for A q 
gives q™ K ~ 2k F u 1 / 2 and A 3 ° x ~ (16/37r)e F u 3 / 2 . These 
approximate locations of g max and A max appear in Fig. [TJi 
with dashed lines, and show a good agreement with the ex- 
act numerical solution of Eq. (j7]) at T = 0. The scaling laws 
correspond to the critical exponents obtained for the equilib- 
rium paramagnetic-ferromagnetic transition with v — 1/2 and 
z = 3 [ 13"]. The temperature dependence of the above quan- 
tities have also been studied numerically and are plotted in 
Fig.[lJ andflj) as dotted lines. We note that the scaling laws 
remain unchanged at finite temperatures. 

Next, we consider the same instability in 2D, which is 
the lower critical dimension. In 2D, there is a quantum 
phase transition to the ferromagnetic state at T — 0, but 
there is no thermodynamic phase transition at finite temper- 
atures. We find, however, that even at nonzero temperatures, 
there still is a critical interaction strength, U C (T), such that 
for Uf > U C (T), there are exponentially growing collective 
modes in the short time dynamics (see the inset of Fig [It). 
At T = 0, Eqs. ^ and Q yield A q D = e F {q/k F )[u/(l + 
m)]a/3 - (q/k F y + 2u + l/(2u+ 1), from which we get 
A™ ax = 2ue F /(2u + 1) ~ 2e P u and cj™ = 
fc Fv /u+ [3+ l/(2u+ l)]/2 ~ V2k F . In Figjlj; and[ljl 
we plot the same quantities as in Fig. [1J and [lp, but for 
the 2D case. The major difference in 2D is that, at finite 
temperatures, q 2 Ja X starts at q = for u = and rapidly 
grows to V2k F as one quenches deeper into the unstable 
side. From a low temperature expansion, we find that at 
T + 0, q^ D ax ~ U 1 /2( T )i/2 exp(ei; ,/2 r ) W hile A max ~ 
M 3 / 2 (Te F ) 1 / 2 exp(e F /2T) for small u. These scaling so- 
lutions go over to the T = scaling forms when u ~ 
(e F /T) exp(— e F /T), which is an exponentially small region 
at small T. 

Pairing instability- In this section we consider the instabil- 
ity of a 3D Fermi gas with attractive interactions at T = 0. 
The case with no spin imbalance has been studied in the liter- 



ature before [14|. Here, we study the general case with spin 
imbalance 031 [T6l . We start by noting that the bare BCS sus- 
ceptibility has a ultraviolet divergence which is regularized by 
expressing the microscopic interaction strength, U /, in terms 
of the s-wave scattering length lTT7l . Following this regular- 
ization procedure, we obtain an expression for the poles of 
the propagator of the pairing fluctuations [the equivalent of 
Eq. 0]: 

m = f d d k / l-n F (&, T )-^-k+ q ,i) m \ 

ATTa s J (2n) d \ ^- a T " C-k +q i + *0+ k 2 J ' 

(8) 

where cj q = J7 q + iA q , and a s < is the final s-wave scat- 
tering length. In Fig. [2^, we plot A q as a function of q, for 
several values of spin imbalance 5 = (ri| — ri|)/(n.| + n^) 
and a final s-wave scattering length k F a s = —0.4. For the 
spin-balanced case we recover the approximate solution A q i=s 
A -2 (? 2 /3fc 2 ,A -C'( (? 4 ), where A = 8e F exp(ir /2a s -2), 
from which we infer g cut = (y / 3/2)Aofc F . We note that 
Aq = Abcs, i-e- the growth rate of the BCS type pairings 
is equal to the BCS gap at equilibrium |[T4l . Thus, for the spin 
balanced case the growth rate decreases monotonically with q 
and the strongest instability is towards pairing with zero center 
of mass momentum. The same picture holds at small but finite 
population imbalances |fl~8). However, as 8 crosses a critical 
value, the maximum of the growth rate continuously shifts to a 
nonzero q. Thus, for large enough 6, the short time dynamics 
is dominated by finite-momentum pairings. This is illustrated 
in Fig.[2]3 where we plot the wave- vector of the most unstable 
state, q max , as a function of spin imbalance S for several val- 
ues of k F a s . As the imbalance is increased further, the zero 
momentum pairing instability disappears completely (corre- 
sponding to Clogston-Chandrasekhar limit |[T9l ). Finally, all 
pairing instabilities disappear at larger population imbalances, 
which corresponds to the sudden end of the curves in Fig.[2j3 
at large <5. 

In the inset of Fig. [2] we have plotted the phase diagram of 
the system as determined by the nature of the dominant un- 
stable mode in the coupling {k F a s ) and the spin imbalance, 
(<5) plane. Note that FFLO is used in a loose sense of finite 
momentum pairing while BCS corresponds to zero momen- 
tum pairing. This diagram resembles the mean-field equilib- 
rium phase diagram [9 |, with the significant difference being 
a larger region dominated by finite-momentum pairings. We 
finally note that instead of directly probing the pairing order 
parameter, one can look for spin correlations associated with 
finite wave-vector pairings through interactions in the Cooper 
channel, which are expected to appear for quenches with suf- 
ficiently large spin imbalances. 

Comment regarding lattice- One interesting extension of 
our work is analysis of the timescale of formation of anti- 
ferromagnetic order in optical lattices for interaction strength 
which are not too large. In these cases, antiferromagnetic or- 
der relies on the BCS type instability of the nested Fermi sur- 
face ll20l . Therefore, the quench dynamics is expected to have 
resemblances with the behavior shown in Fig. [2] 
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FIG. 2: Unstable collective modes in the BCS instability in 3D. (a) 
Growth rate of modes, A q , as a function of wavevector q for vari- 
ous values of the population imbalance 8 = (top curve), 0.0025, 
. . . , 0.0225 (bottom curve); (fc_pa s = —0.4). (b) The most unstable 
wavevector g max as a function of S for various values of the kpa a (in- 
dicated). As & increases, g milx goes through a continuous transition; 
the curves terminate at large 5, above which all pairing modes are 
stable. Inset: the phase diagram in the kpa s -8 plane, showing the 
character of the most unstable mode. 



Cold atom experiments- Finally, we would like to comment 
on the possibility of observing the described dynamics of tex- 
ture growth in realistic cold atom experiments. There are two 
issues here: (a) the length scale corresponding to the most un- 
stable mode and (b) the time scale associated to the growth of 
unstable modes. The length scale l/qWx must be sufficiently 
large so that patterns can be resolved. Achieving large pattern 
sizes requires a quench to a final coupling which is close to the 
critical value. However, the growth rate A max is small in this 
regime and the textures develop over a long time. Therefore, 
one needs to find a window of Uf, where the length scale is 
resolved in experiments while the timescale is short enough to 
prevent considerable atom losses. In recent experiments Q, 
the quench time was ~ 4 ms, which was followed by a hold 
time of upto ~ 12 ms; the imaging resolution was ~ lOfc^ 1 . 

We note that the quench approximation is applicable for 
times shorter than the corresponding (time dependent) inverse 
growth rate, i.e. A max [t/(i)] _1 ~ t— t c , where t c is the time at 
which the instability is crossed, U(t c ) = U c . At longer times, 
the evolution of the magnetization is dominated by the much 
slower dynamics of the coalescence of the magnetic domains. 
Assuming that U (t) in the MIT experiments is ramped up in 
time from to Ut at a constant rate, the quench approxima- 
tion is applicable to the MIT experiment with final couplings 
of Uf/U c < 1.06. For data points where U f /U c > 1.06, the 
quench dynamics stops at the point U(t)/U c ~ 1.06 which 
results in a pattern of magnetic domains with typical size 
~ 2k p 1 that coarsens very slowly over the remainder of the 
experment. In the MIT data, the typical step size in the final 
interaction strength was Akpa s ~ 0.3, thus we would ex- 
pect that the smallest U f that exceeded U c was Uf ~ 1.06 U c . 
Therefore, we would not expect to see a texture with a wave- 
length longer than ~ 2/c^ 1 . This is consistent with the fact 
that individual textures were not resolved in the current set of 
experiments, as the MIT setup can only resolve textures with 



wavelength of ~ 10k p 1 . However, our theory indicates that 
domains of this size could be achieved if the quench is stopped 
at U f/U c ~ 1.004, after a hold time of ~ 14ms, which could 
be possible with current experiments via a finer tuning of the 
final interaction strength (although the inhomogeneous den- 
sity profile in the trap and three body recombinations do make 
it difficult to observe). 
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